SCP - Centrality
1 Carga de librerías
if (!(require(DT)))
install.packages("DT")
library(DT)
if (!require(readxl))
install.packages("readxl")
library(readxl)
# Rfast para matriz de varianza combinada
if (!requireNamespace("Rfast"))
install.packages("Rfast")
library(Rfast)
if (!requireNamespace("ggplot2"))
install.packages("ggplot2")
library(ggplot2)
if (!requireNamespace("plotly"))
install.packages("plotly")
library(plotly)
# GridFCM
library(devtools)
if (!requireNamespace("GridFCM.practicum"))
install_github("asanfe/GridFCM.practicum", quietly = TRUE)
library(GridFCM.practicum)
# Viridislilte
if (!requireNamespace("viridisLite"))
install.packages("viridisLite")
library(viridisLite)
# Test para normalidad multivariante
if (!requireNamespace("MVN"))
install.packages("MVN")
library(MVN)
if (!requireNamespace("ggpattern", quietly = TRUE))
install.packages("ggpattern")
library(ggpattern)
if (!requireNamespace("factoextra", quietly = TRUE))
install.packages("factoextra")
library(factoextra)
if (!requireNamespace("cluster", quietly = TRUE))
install.packages("cluster")
library(cluster)
if (!requireNamespace("RColorBrewer", quietly = TRUE))
install.packages("RColorBrewer")
library(RColorBrewer)
if (!requireNamespace("rcartocolor", quietly = TRUE))
install.packages("rcartocolor")
library(rcartocolor)
if (!requireNamespace("dplyr", quietly = TRUE))
install.packages("dplyr")
library(dplyr)2 Importación y resumen
2.1 Importación del objeto RDA
# Objetos de sesión de ejemplo de la PEC
path <- '../CentralityTest/data.RData'
load(path)
samples.raw.df <- data
samples.df <- data.frame(ID = samples.raw.df$dataset$ID,
gender = samples.raw.df$dataset$gender,
age = samples.raw.df$dataset$age,
edu = samples.raw.df$dataset$edu,
status = "Error")
for(i in 1:nrow(samples.df)) {
id <- samples.df$ID[i]
tryCatch({
wimp <- data$grids[[id]]$WT # Accede al wimp asociado al ID
wphm <- GridFCM.practicum::ph_index(wimp = wimp, method = "wnorm", std = "none")
samples.df$status[i] <- "Ok"
}, error = function(e){
cat("Error procesando ID:", id, "\n")
})
}## Error procesando ID: P0166567
## Error procesando ID: P0512115
## Error procesando ID: P0606903
## Error procesando ID: P0666512
## Error procesando ID: P0870223
## Error procesando ID: P0910424
## Error procesando ID: P1102149
## Error procesando ID: P1123165
## Error procesando ID: P1140623
## Error procesando ID: P1312902
## Error procesando ID: P1426704
## Error procesando ID: P1554581
## Error procesando ID: P1891931
## Error procesando ID: P1933446
# Calcula las frecuencias de status de procesamiento
freqs.status <- table(samples.df$status, useNA = "no")
# Calcula los porcentajes
percent.status <- prop.table(freqs.status) * 100
results.summary <- data.frame(
ResultadoCarga = names(freqs.status),
Casos = as.integer(freqs.status),
Porcentaje = round(100 * prop.table(freqs.status), 3)
)
results.summary <- results.summary[, c("ResultadoCarga", "Casos", "Porcentaje.Freq")]2.2 Resultado de procesamiento
2.4 Wimp del sujeto a presentar
#path <- 'Wimp_Ejemplo.xlsx'
#opr <- TRUE
#wimp <- GridFCM.practicum::importwimp(path = path, opr = opr, sheet = 1)
# Objetos de sesión de ejemplo de la PEC
path <- '../CentralityTest/data.RData'
load(path)
samples.raw.df <- data
# Selección de un sujeto
id.sample <- params$id.sujeto.entrada
wimp <- data$grids[[id.sample]]$WT
sujeto.df <- data$dataset[data$dataset$ID == id.sample,]
bertin(wimp$openrepgrid, colors = c("palegreen", "darkgreen"))3 Exploración de la Wimp
3.3 E/S de los constructos. Método Simple
3.4 E/S de los constructos. Método wnorm
4 Distancia de Mahalanobis y distribución de datos
4.1 Test de Mardia para análisis multivariante
Llevamos a cabo previamente un test de Mardia para constrastar la normalidad multivariante de los datos, a fin de determinar la pertinencia del punto de corte basado en adecuación a distribución Chi-cuadrado de distancia de Mahalanobis
## $multivariateNormality
## Test Statistic p value Result
## 1 Mardia Skewness 8.4418086965317 0.0766708484984564 YES
## 2 Mardia Kurtosis 0.152819936911057 0.878540275026683 YES
## 3 MVN <NA> <NA> YES
##
## $univariateNormality
## Test Variable Statistic p value Normality
## 1 Anderson-Darling p 0.3197 0.4706 YES
## 2 Anderson-Darling h 0.5754 0.0998 YES
##
## $Descriptives
## n Mean Std.Dev Median Min Max 25th
## p 10 6.268880e-01 0.07552411 0.60218060 0.5110669 0.7507256 0.5802010
## h 10 -1.109681e-17 0.21125846 -0.07780234 -0.2527941 0.4152149 -0.1307287
## 75th Skew Kurtosis
## p 0.6672639 0.2886611 -1.248914
## h 0.1341522 0.7127256 -1.001594
4.2 Test de resultado de la función
4.3 Distribución Chi-cuadrado
# Definimos los grados de libertad para la distribución chi-cuadrado
df <- 2
# Generamos los valores de la distribución
x <- seq(qchisq(0.001, df), qchisq(0.999, df), length.out = 1000)
y <- dchisq(x, df)
# Calculamos los puntos de corte para el 20% superior
sign.level <- 0.2
cut_high <- qchisq(1- sign.level, df)
# Dataframe para la gráfica
data <- data.frame(x = x, y = y)
ggplot(data, aes(x = x, y = y)) +
geom_line() +
geom_ribbon(data = data %>% filter(x > cut_high), aes(ymax = y), ymin = 0, fill = 'salmon', alpha = 0.5) +
geom_vline(xintercept = cut_high, color = "red", linetype = "dashed") +
labs(title = 'Distribución Chi-cuadrado con puntos de corte del 80%', x = 'Valor', y = 'Densidad') +
theme_minimal()4.4 Gráfica de barras de distancias de Mahalanobis y punto de corte
# Distancia de Mahalanobis
test.bp.wmahalanobis <- GridFCM.practicum::mahalanobis_index(wimp = wimp, method = "wnorm", std = FALSE)
test.wmahalanobis.df <- as.data.frame(test.bp.wmahalanobis)
# Colores de los constructos
#test.wmahalanobis.df$constructo <- rownames(test.wmahalanobis)
test.wmahalanobis.df$constructo <- wimp$constructs$left.poles
# Valoración del ideal
test.wmahalanobis.df$idealdirect <- wimp$ideal$direct
# Columna para identificar constructos dilemáticos
#test.wmahalanobis.df$fill.color <- ifelse(test.wmahalanobis.df$idealdirect == 4, "yellow2", "honeydew")
test.wmahalanobis.df$fill.color <- construct_colors(wimp= wimp, mode = "red/green")
# Ordenamos las barras en orden decreciente
test.wmahalanobis.df <- test.wmahalanobis.df %>%
arrange(desc(m.dist))
# Convertimos 'constructo' en un factor con los niveles en el orden deseado
test.wmahalanobis.df$constructo <- factor(test.wmahalanobis.df$constructo, levels = test.wmahalanobis.df$constructo)
# Punto de corte distribución Chi-Cuadrado
sign.level <- 0.2
df <- ncol(test.wphm)
chi.square.cutoff <- qchisq(1 - sign.level, df)
#media_m_dist <- mean(test.wmahalanobis.df$m.dist)4.5 Constructos supraordenados
# Crear el histograma de constructos supraordenados
# Filtramos por los constructos donde el valor de 'h' es mayor que cero
test.wmahalanobis.df.sup <- test.wmahalanobis.df %>%
filter(h > 0)
bar_plot <- ggplot(test.wmahalanobis.df.sup, aes(x = constructo, y = m.dist, fill = fill.color)) +
geom_bar(stat = "identity", color = "black", linewidth = 0.25) +
geom_hline(yintercept = chi.square.cutoff, linetype = "dashed", color = "darkgreen", linewidth = 1) +
scale_fill_identity() + # Usa los colores asignados directamente
theme_minimal() +
theme(
axis.text.x = element_text(angle = 45, hjust = 1),
legend.position = "none" # Ocultar leyenda para fill
) +
labs(x = "Constructos con h > 0", y = "Distancia de Mahalanobis", title = "Constructos con h>0 por distancia de Mahalanobis")
# Mostramos el gráfico
print(bar_plot)4.6 Constructos subordinados
# Crear el histograma de constructos subordinados
# Filtramos por los constructos donde el valor de 'h' es menor que cero
test.wmahalanobis.df.sub <- test.wmahalanobis.df %>%
filter(h < 0)
bar_plot <- ggplot(test.wmahalanobis.df.sub, aes(x = constructo, y = m.dist, fill = fill.color)) +
geom_bar(stat = "identity", color = "black", linewidth = 0.25) +
geom_hline(yintercept = chi.square.cutoff, linetype = "dashed", color = "darkgreen", linewidth = 1) +
scale_fill_identity() + # Usa los colores asignados directamente
theme_minimal() +
theme(
axis.text.x = element_text(angle = 45, hjust = 1),
legend.position = "none" # Ocultar leyenda para fill
) +
labs(x = "Constructos con h < 0", y = "Distancia de Mahalanobis", title = "Constructos con h<0 por distancia de Mahalanobis")
# Mostramos el gráfico
print(bar_plot)4.7 Distribución de valores en P
4.7.1 Distribución de valores de P
# Crear el histograma de valores en P
test.wmahalanobis.df.sortP <- test.wmahalanobis.df %>%
arrange(desc(p))
mean.p <- mean(abs(test.wmahalanobis.df.sortP$p))
# Convertir 'constructo' en un factor para mantener el orden en el gráfico
test.wmahalanobis.df.sortP$constructo <- factor(test.wmahalanobis.df.sortP$constructo,
levels = test.wmahalanobis.df.sortP$constructo)
bar_plot <- ggplot(test.wmahalanobis.df.sortP, aes(x = constructo, y = p, fill = fill.color)) +
geom_bar(stat = "identity", color = "black", linewidth = 0.25) +
geom_hline(yintercept = mean.p, linetype = "dashed", color = "darkgreen", linewidth = 1) +
scale_fill_identity() + # Usa los colores asignados directamente
theme_minimal() +
theme(
axis.text.x = element_text(angle = 45, hjust = 1),
legend.position = "none" # Ocultar leyenda para fill
) +
labs(x = "Constructos", y = "Valor de P", title = "Constructos por distancia en P")
# Mostramos el gráfico
print(bar_plot)4.7.2 Test Shapiro-Wilk (muestras pequeñas o moderadas) para variable P
# Test de normalidad de Saphiro-Wilk
norm.test <- shapiro.test(test.wmahalanobis.df$p)
# Imprime el resultado
print(norm.test)##
## Shapiro-Wilk normality test
##
## data: test.wmahalanobis.df$p
## W = 0.9447, p-value = 0.6064
De acuerdo con el resultado de la prueba, la normalidad de la distribución de los valores de P es: TRUE
4.7.3 Gráfica Cuantil-Cuantil
datos <- test.wmahalanobis.df$p
# Gráfica q-q para comprobar la normalidad
qqnorm(datos)
qqline(datos, col = "red")4.7.4 Punto de corte basado en distribución normal de P
# Media y desviación típica de la distribución
mean.p <- mean(test.wmahalanobis.df$p)
sd.p <- sd(test.wmahalanobis.df$p)
# Definir el rango de valores para X basado en la media y desviación típica
x.values <- seq(from = mean.p - 4 * sd.p, to = mean.p + 4 * sd.p, length.out = 1000)
# Crear un dataframe con los valores de X y la densidad de una distribución normal para esos valores
norm.df <- data.frame(x = x.values, y = dnorm(x.values, mean = mean.p, sd = sd.p))
# Punto de corte
cut.low <- qnorm(0.15, mean = mean.p, sd = sd.p)
plot <- ggplot(norm.df, aes(x = x, y = y)) +
geom_line() +
geom_vline(xintercept = cut.low, linetype = "dashed", color = "red", linewidth = 1) +
geom_vline(xintercept = mean.p, linetype = "dashed", color = "blue", linewidth = 1) +
geom_area(data = subset(norm.df, x <= cut.low), fill = "lightblue", alpha = 0.2) +
theme_bw() +
theme(
panel.grid.major = element_line(linewidth = 0.5, linetype = 'solid', colour = "lightgrey"),
panel.grid.minor = element_blank(),
legend.position = "none"
) +
scale_x_continuous(name = "Valor de P") +
scale_y_continuous(name = "Densidad") +
labs(title = paste("Distribución Normal con Media en", round(mean.p, 2),
"y Punto de Corte en", round(cut.low, 2)))
# Mostrar la gráfica
print(plot)5 Centralidad y orden subjetivo de los constructos
# Crear una nueva columna para almacenar el valor de importancia en test.wmahalanobis.df
test.wmahalanobis.df$importanciaSubjetiva <- NA
test.wmahalanobis.df$totalCentrl <- NA
# Asignar el valor de cen.ord.cX.test basado en el número de fila
for (i in 1:nrow(test.wmahalanobis.df)) {
constructo.actual <- test.wmahalanobis.df$constructo[i]
# Buscar este constructo en 'sujeto.df' desde c1l a c10l
for (j in 1:10) {
col.constructo <- paste("c", j, "l", sep = "")
col.importanciaSub <- paste("cen.ord.c", j, ".test", sep = "")
col.totalCentrl <- paste("cen.c", j, ".test", sep = "")
if (constructo.actual %in% sujeto.df[[col.constructo]]) {
# Si el constructo se encuentra, asignar la importancia subjetiva correspondiente
test.wmahalanobis.df$importanciaSubjetiva[i] <- sujeto.df[[col.importanciaSub]]
test.wmahalanobis.df$totalCentrl[i] <- sujeto.df[[col.totalCentrl]]
break # Salir del bucle interno una vez que se asigna el valor
}
}
}
DT::datatable(select(test.wmahalanobis.df, p, h, m.dist, hub, importanciaSubjetiva, totalCentrl))6 Representación en espacio P-H
6.1 Sin estandarización - plotly sin marcar área no viable ni constructos centrales
6.2 Espacio PH con coloreado de área no útil y marcado de outliers. Función de representación
6.2.1 Sin estandarización
6.2.2 Con estandarización basada en aristas
6.2.3 Sin estandarización, marcando área de coordenadas no viables. Sin puntos
wp1.grph <- GridFCM.practicum::graph_ph(wimp = wimp, method = "wnorm", std = "none", sign.level = 0.2,
mark.nva = TRUE, mark.hub = TRUE, show.points = FALSE)
wp1.grph## No trace type specified:
## Based on info supplied, a 'scatter' trace seems appropriate.
## Read more about this trace type -> https://plotly.com/r/reference/#scatter
## No scatter mode specifed:
## Setting the mode to markers
## Read more about this attribute -> https://plotly.com/r/reference/#scatter-mode
7 Otros métodos de centralidad
7.1 Eigenvectores sobre matriz de implicaciones
La puntuación de centralidad para cada constructo se obtiene sumando los cuadrados de las componentes del primer y segundo vector propio, ponderados por los respectivos valores propios. Este vector de puntuciones representa cuánto contribuye cada constructo a las dos principales dimensiones de variabilidad en los datos.
eigen_index <- function(wimp) {
# Matriz de adyacencia
adj.matrix <- wimp$scores$implications
# Vectores y valores propios
results <- eigen(adj.matrix)
# Emplearemos información de los dos primeros vectores propios
eigenvectorF <- results$vectors[,1]
eigenvectorS <- results$vectors[,2]
# Puntuación combinada
centralidad <- eigenvectorF^2 * results$values[1] + eigenvectorS^2 * results$values[2]
# Dataframe para resultados
df.centrality <- data.frame(
constructs = wimp$constructs$constructs,
leftpoles = wimp$constructs$left.poles,
rightpoles = wimp$constructs$right.poles,
centrality = abs(round(centralidad, 3))
)
return(df.centrality)
}Resultado de la función:
7.2 Análisis de componentes principales (PCA)
La siguiente función calcula la centralidad de los constructos utilizando un PCA aplicado a la matriz de implicaciones. Se calcula la varianza explicada por cada componente principal y se usa para ponderar las cargas de los dos primeros componentes. La centralidad de cada constructo se determina sumando las cargas cuadradas de los dos primeros componentes principales, ponderadas por la varianza explicada por estos componentes.
pca_index <- function(wimp) {
# Foco del PCA
adj.matrix <- wimp$scores$implications
# Análisis de componentes principales, y varianza explicada por cada componente
pca.result <- prcomp(adj.matrix, center = TRUE, scale = TRUE)
explained.variance <- pca.result$sdev^2 / sum(pca.result$sdev^2)
# Cargas de los componentes
#loadings.pca <- pca.result$rotation
# Suma de cargas en los dos primeros componentes principales
loadings.add <- rowSums((pca.result$rotation[, 1:2]^2) * explained.variance[1:2])
# Creamos dataframe con los nombres de los constructos y las puntuaciones de la suma de cargas en PC1 a PC2
pca.df <- data.frame(
constructs = wimp$constructs$constructs,
leftpoles = wimp$constructs$left.poles,
rightpoles = wimp$constructs$right.poles,
centrality = abs(round(loadings.add, 3))
)
return(pca.df)
}Resultado de la función:
7.3 PCA Plot
# Foco del PCA
adj.matrix <- wimp$scores$implications
# Análisis de componentes principales
pca.result <- prcomp(adj.matrix, center = TRUE, scale = TRUE)
# Extraer los dos primeros componentes principales
pca.comp <- as.data.frame(pca.result$x[, 1:2])
pca.comp$constructs <- wimp$constructs$constructs
# Crea la gráfica de dispersión
pca.plot <- plot_ly(data = pca.comp, x = ~PC1, y = ~PC2, type = 'scatter', mode = 'markers',
hoverinfo = 'text+x+y',
marker = list(size = 10, opacity = 0.8)) %>%
layout(title = 'PCA de matriz de adyacencia',
xaxis = list(title = 'PCA 1'),
yaxis = list(title = 'PCA 2'),
hovermode = 'closest',
plot_bgcolor = "white",
font = list(family = "Arial"),
showlegend = FALSE) %>%
# Add annotations for each point
add_annotations(data = pca.comp, x = ~PC1, y = ~PC2, text = ~constructs,
showarrow = FALSE, xanchor = 'center', yanchor = 'bottom', font = list(size = 12))
# Muestra la gráfica
pca.plot7.4 Comparativa de métodos de centralidad
En este punto, ofrecemos un “ranking” de centralidad de los constructos según los tres métodos empleados.
# Ejecución de las fuciones sobre la wimp del caso
mahalanobis.res.lst <- mahalanobis_index(wimp)
eigen.res <- eigen_index(wimp)
pca.res <- pca_index(wimp)
mahalanobis.res <- as.data.frame(mahalanobis.res.lst)
mahalanobis.res$constructs <- rownames(mahalanobis.res)
unq.constructs <- unique(c(mahalanobis.res$constructs, eigen.res$constructs, pca.res$constructs))
# Preparar un data frame final
comb.df <- data.frame(constructs = unq.constructs)
comb.df <- merge(comb.df, mahalanobis.res[, c("constructs", "m.dist")], by = "constructs", all.x = TRUE)
comb.df <- merge(comb.df, eigen.res[, c("constructs", "centrality")], by = "constructs", all.x = TRUE,
suffixes = c(".mahalanobis", ".eigen"))
comb.df <- merge(comb.df, pca.res[, c("constructs", "centrality")], by = "constructs", all.x = TRUE)
# Redondea las columnas numéricas a 3 decimales
comb.df <- comb.df %>%
mutate(across(where(is.numeric), round, digits = 3))## Warning: There was 1 warning in `mutate()`.
## ℹ In argument: `across(where(is.numeric), round, digits = 3)`.
## Caused by warning:
## ! The `...` argument of `across()` is deprecated as of dplyr 1.1.0.
## Supply arguments directly to `.fns` through an anonymous function instead.
##
## # Previously
## across(a:b, mean, na.rm = TRUE)
##
## # Now
## across(a:b, \(x) mean(x, na.rm = TRUE))
# Renombramos columnas
names(comb.df)[2:4] <- c("Centr_Mahalanobis", "Centr_Eigen", "Centr_PCA")
# Ranking de los constructos dentro de cada método
comb.df$Rank_Mahalanobis <- rank(-comb.df$Centr_Mahalanobis, ties.method = "first")
comb.df$Rank_Eigen <- rank(-comb.df$Centr_Eigen, ties.method = "first")
comb.df$Rank_PCA <- rank(-comb.df$Centr_PCA, ties.method = "first")
comb.df <- comb.df[order(comb.df$Rank_Mahalanobis),]
DT::datatable(comb.df)7.5 Combinación con medidas de importancia subjetiva
mahalanobis.res$importanciaSubjetiva <- NA
mahalanobis.res$totalCentrl <- NA
for (i in 1:nrow(mahalanobis.res)) {
mahalanobis.res$importanciaSubjetiva[i] <- sujeto.df[[paste("cen.ord.c", i, ".test", sep = "")]]
mahalanobis.res$totalCentrl[i] <- sujeto.df[[paste("cen.c", i, ".test", sep = "")]]
}
# Dataframe combinado para resultados
centrality.comb.df <- merge(mahalanobis.res, comb.df, by = "constructs", all = TRUE)
# Redondea las columnas numéricas a 3 decimales
centrality.comb.pres.df <- centrality.comb.df %>%
mutate(across(where(is.numeric), round, digits = 3))
DT::datatable(centrality.comb.pres.df)8 Análisis por conglomerados
8.1 Determinación de número óptimo de conglomerados
El cálculo se basará en la distancia de Mahalanobis entre pares de puntos. La distancia de Mahalanobis se calcula utilizando la fórmula:
\[d(\mathbf{x}, \mathbf{y}) = \sqrt{(\mathbf{x} - \mathbf{y})^T \mathbf{S}^{-1} (\mathbf{x} - \mathbf{y})}\]
Donde:
- \(\mathbf{x}\) y \(\mathbf{y}\) son los dos vectores que representan los dos puntos en el espacio.
- \(\mathbf{S}^{-1}\) es la matriz inversa de la matriz de covarianza de los datos.
## [1] 4
Tenemos un número máximo de 4 conglomerados en nuestros datos.
8.2 Representación de números de conglomerados óptimo
Adecuación de cohesión y separación de cada punto según pertenezca a distintos conglomerados:
# Lista que albergará las distintas gráficas de silueta
lista.graf.sil <- list()
# Valores intermedios que ya calcula .optimal.num.clusters
max.clusters <- length(wimp$constructs$constructs) - 1
ph.mat <- GridFCM.practicum::ph_index(wimp = wimp, method = "wnorm", std = "none")
#rownames(test.dist) <- as.character(wimp$constructs$right.poles)
#distancias <- dist(test.dist, method = "euclidean")
#------------------------------
# Matriz de disimilaridad modelada como matriz de distancias de Mahalanobis
# Matriz de covarianzas
cov.matrix <- cov(ph.mat) # Calcula la matriz de covarianza
# Vector de medias
means.vector <- colMeans(ph.mat)
# Inicializa una matriz para guardar las distancias de Mahalanobis
n <- nrow(ph.mat)
dist.mat <- matrix(NA, n, n)
# Calcula la distancia de Mahalanobis entre cada par de filas en ph.mat
for (i in 1:n) {
for (j in i:n) {
diff <- ph.mat[i, ] - ph.mat[j, ]
dist.mat[i, j] <- sqrt(t(diff) %*% solve(cov.matrix) %*% diff)
dist.mat[j, i] <- dist.mat[i, j] # La matriz es simétrica
}
}
# Hacemos 0 en la diagonal
diag(dist.mat) <- 0
row.names(dist.mat) <- row.names(ph.mat)
colnames(dist.mat) <- row.names(ph.mat)
#---------------------------
# Preparamos una lista con diversas representaciones gráficas de siluetas (de 2 a 10 clústeres)
for(j in 2:min(13, max.clusters)){
it.pam <- cluster::pam(dist.mat, j, diss = TRUE)
p <- factoextra::fviz_silhouette(it.pam, label = FALSE, print.summary = FALSE)
lista.graf.sil[[j-1]] <- p
}
# Organizar los gráficos en una matriz de 4x3, y los presentamos
gridExtra::grid.arrange(grobs = lista.graf.sil, ncol = 3, nrow = 4)8.2.1 Dendrograma
## Registered S3 method overwritten by 'dendextend':
## method from
## text.pvclust pvclust
8.2.2 ClusPlot
act.cex <- par("cex")
par(cex = 0.8)
# Calculamos el objeto de partición PAM para los k clústeres obtenidos anteriormente
opt.pam <- cluster::pam(dist.mat, k, diss = TRUE)
# Colores de los clusters
clus.colors <- carto_pal(n = k, "Peach")
cluster::clusplot(x = dist.mat,
clus = opt.pam$clustering,
shade = TRUE,
color = TRUE,
col.clus = clus.colors,
col.p = 'darkblue',
diss = TRUE,
labels = 3)